The invention relates to the field of attenuating noise in signals. More particularly, the invention relates to attenuating incoherent, impulsive and coherent noise in seismic data. Even more particularly, the invention relates to using a weighted slant stack to attenuate incoherent, impulsive and coherent noise in seismic data.
Geophysicists collect seismic data and analyze it to determine the characteristics of underground and undersea formations. Such information is useful in the search for hydrocarbons and other minerals.
In its raw, unprocessed form, seismic data consists of a large number, typically millions, of xe2x80x9ctraces.xe2x80x9d Each trace is a recording of a signal which represents the seismic energy emitted by a seismic transmitter, reflected off underground formations and received by a receiver. The seismic energy may be transmitted by an underground or undersea explosion or the action of a vibrator truck or through some other means of imparting energy into the earth. The seismic data typically includes signals gathered from receivers arrayed in two dimensions over a geographical region to be analyzed. Typically, a number of traces is gathered from each receiver, with the location of the source of the seismic energy being varied from one trace to the next.
In many cases, analysis of the seismic data is complicated by unwanted noise recorded in the traces along with the signal. The noise may include incoherent noise, impulsive noise or coherent noise or all of those together. Incoherent noise is random noise such that the incoherent noise from two or more traces tends to have a low level of correlation. Impulsive noise tends to appear in the form of spikes or abnormally large amplitudes in isolated portions of a seismic trace or in a whole trace in a number of traces of the seismic data set. Coherent noise tends to correlate among traces and, in many cases, can be grouped early in the processing of the seismic data into sets having identifiable xe2x80x9ctime dips.xe2x80x9d
Incoherent noise is traditionally attenuated by xe2x80x9cstackingxe2x80x9d or adding a group of traces together. Since incoherent noise is random, it tends to be eliminated by the stacking process. More effective elimination of incoherent noise is accomplished by a technique called f-x-y prediction error filtering (fxy PEF). This technique predicts the signal in the x-y space at each temporal frequency f, thereby reducing the incoherent noise, which is not predictable.
Impulsive noise is usually attenuated using statistical techniques in running windows in which abnormally large amplitudes are edited based on certain types of average measurement. Other techniques use adaptive filtering or neural networks. These algorithms are based on prior training on data containing a typical signal and typical noise. Once the algorithm understands what the signal is and what the noise is, the algorithm is applied to a data set.
Existing techniques for attenuating coherent noise include transforming the traces from the space-time domain into a transform domain, such as the f-k or Radon or slant-stack domain, in which the signal is separated from the noise. The noise is muted in the transform domain and the signal is transformed back into the space-time domain for further processing. Since the noise was muted in the transform domain, it is attenuated in the space-time domain. For example, coherent noise having a particular time dip may be separated from the signal in transforming a trace into the Radon or slant-stack transform domain where it can be muted.
Other existing technologies include combining the diversity stack technique with the slant stack technique. In this technique, incoherent noise is attenuated by a variation on the diversity stack technique and coherent noise is attenuated by transforming the traces to the Radon or slant-stack domain, muting the traces associated with the coherent noise, and transforming the traces back into the space-time domain. The impulsive noise is also attenuated by the diversity technique by suppressing abnormally large amplitudes present in the amplitude array to be stacked.
In general, in one aspect, the invention features a method of attenuating noise in seismic data including a plurality of input traces. The method includes transforming the seismic data from the space-time domain into the slant-stack domain. Seismic data having a preselected characteristic is excluded when the transforming into the slant-stack domain. The transformed data is inverse transformed from the slant-stack domain into the time space domain.
Implementations of the invention may include one or more of the following. Transforming may include transforming the seismic data from the space-time domain into the xcfx84xe2x88x92p domain. Excluding may include excluding data corresponding to the coordinates: t=xcfx84+px, where t is time, xcfx84 is a t-axis intercept, x is a signed distance from an origin and p is a slope. Transforming may include transforming the seismic data from the space-time domain into the xcfx84xe2x88x92pxxe2x88x92py domain. Excluding may include excluding data corresponding to the coordinates: t=xcfx84+pxx+pyy, where t is time, xcfx84 is a t-axis intercept, x is a signed distance from an origin in a first direction, y is a signed distance from an origin in a second direction, px is slope in the first direction, and py is slope in the second direction.
In general, in another aspect, the invention features a method of attenuating noise in seismic data. The method includes accepting seismic data comprising a plurality of traces arrayed in two physical dimensions and a time dimension. The method further includes transforming the seismic data into the slant stack domain using a Radon transform.
Implementations of the invention may include one or more of the following. The Radon transform may include a three-dimensional transform using the following equation:       T    ⁡          (              τ        ,                  p          x                ,                  p          y                    )        =                                          N            ·                                          ∑                                  x                  =                  a                                b                            ⁢                                                ∑                                      y                    =                    c                                    d                                ⁢                                                                            S                                              x                        ,                        y                                                              ⁡                                          (                                                                        τ                          +                                                                                    p                              x                                                        ⁢                            x                                                    +                                                                                    p                              y                                                        ⁢                            y                                                                          ,                        x                        ,                        y                                            )                                                        ·                                                                                                                    F                              x                ,                y                                      ⁡                          (                                                τ                  +                                                            p                      x                                        ⁢                    x                                    +                                                            p                      y                                        ⁢                    y                                                  ,                x                ,                y                            )                                                          ∑                  x          =          a                b            ⁢                        ∑                      y            =            c                    d                ⁢                              F                          x              ,              y                                ⁡                      (                                          τ                +                                                      p                    x                                    ⁢                  x                                +                                                      p                    y                                    ⁢                  y                                            ,              x              ,              y                        )                              
where:
N is the number of traces in the plurality of traces;
Sx,y(t) is a subset of the plurality of traces;
Fx,y(t) is a scaling function;
x is a distance in a first direction;
px is a slope in the first direction;
y is a distance in a second direction;
py is a slope in the second direction;
a and b for the first direction and c and d for the second direction define a volume to be transformed;
t is time in the space-time domain;
xcfx84 is intercept time in the xcfx84xe2x88x92p domain; and
T(xcfx84,px, py) is the output data.
The following equation may apply:             F              x        ,        y              ⁢          (      t      )        =      1                            S          _                          x          ,          y                    ⁢              (        t        )            
where
{overscore (S)}x,y(t) is an average of Sx,y(t) over predetermined ranges of x, y and t for each x, y and t.
The following equation may apply:             F              x        ,        y              ⁢          (      t      )        =      1                            S          _                          x          ,          y                2            ⁢              (        t        )            
where
{overscore (S)}2x,y(t) is an average of S2x,y(t) over predetermined ranges of x, y and t for each x, y and
a, b, c and d may be such that x and y span Sx,y(t), a, b, c and d may be such that x and y define a subset of Sx,y(t). The subset defining Sx,y(t) may be the full set.
The method may include Fourier transforming T(xcfx84,px, py), multiplying the Fourier-transformed data by xcfx892; and inverse-Fourier transforming the Fourier-transformed data.
The Radon transform may include a two-dimensional transform using the following equation:       T    ⁢          (              τ        ,        p            )        =            N      ·                        ∑                      x            =            a                    b                ⁢                                            S                              x                ,                y                                      ⁢                          (                                                τ                  +                                      p                    ⁢                                          xe2x80x83                                        ⁢                    x                                                  ,                x                            )                                ·                                    F                              x                ,                y                                      ⁢                          (                                                τ                  +                                      xe2x80x83                                    ⁢                                      p                    ⁢                                          xe2x80x83                                        ⁢                    x                                                  ,                x                            )                                                          ∑                  x          =          a                b            ⁢                        F                      x            ,            y                          ⁢                  (                                    τ              +                              p                ⁢                                  xe2x80x83                                ⁢                x                                      ,            x                    )                    
where:
N is the number of traces in the plurality of traces;
Sx,y(t) is a subset of the plurality of traces;
Fx,y(t) is a scaling function;
x is a distance from the origin;
p is a slope;
a and b define a slice of data to be transformed;
t is time in the space-time domain;
xcfx84 is intercept time in the xcfx84xe2x88x92p domain; and
T(xcfx84,p) is the output data.
The method may include Fourier transforming T(xcfx84,p), multiplying the Fourier-transformed data by xcfx89, and inverse-Fourier transforming the Fourier-transformed data.
The method may include applying an anti-aliasing filter. Applying may be performed in the slant stack domain.
The method may include applying a matched filter to compensate for limited spatial aperture.
The method may include inverse Radon transforming the transformed data.
In general, in another aspect, the invention features a method of attenuating noise in seismic data. The seismic data includes a plurality of traces arrayed in two physical dimensions and the time dimension. The method includes:
(a) selecting a volume, the volume comprising a subset of the plurality of traces;
(b) transforming the volume into the slant stack domain using a Radon transform;
(c) inverse transforming the transformed volume into the time domain using a inverse Radon transform; and
(d) repeating elements (b) and (c) with a new volume.
Implementations of the invention may include one or more of the following. The new volume may overlap the volume.
The method may include iteratively repeating elements (b) and (c) with new volumes for each iteration, each successive new volume overlapping the preceding new volume.
The plurality of traces may be arrayed in rows and the overlap between each successive new volume and the preceding new volume may be one row.
In general, in another aspect, the invention features a method of attenuating noise in zero offset and/or common offset seismic data. The seismic data includes a plurality of input traces. The method includes:
(a) computing scaling traces from the input traces;
(b) modifying the input traces using the scaling traces;
(c) slant-stack transforming the modified input traces, excluding slopes corresponding to coherent noise from the slant-stack transform, to produce slant-stack-domain signal traces;
(d) slant-stack transforming the scaling traces to produce slant-stack-domain scaling traces;
(e) modifying the slant-stack-domain signal traces using the slant-stack-domain scaling traces;
(f) scaling and filtering the modified slant-stack-domain signal traces to produce weighted slant stacked traces;
(g) inverse slant-stack transforming the weighted slant stacked traces; and
(h) scaling the inverse-slant-stack-transformed weighted slant stacked traces.
In general, in another aspect, the invention features a computer program for attenuating noise in seismic data. The computer program executes on a computer processor. The computer program resides on computer readable media. The computer program includes computer code for causing the computer processor to transform the seismic data from the time-space domain into the slant-stack domain. The computer code for causing the computer processor to transform the seismic data includes computer code for excluding seismic data having a preselected characteristic and transforming the transformed data from the slant-stack domain into the time-space domain.
In general, in another aspect, the invention features an apparatus for attenuating noise in seismic data. The seismic data includes a plurality of input traces. The apparatus includes a means for transforming the seismic data from the time-space domain into the slant-stack domain. The transforming means includes means for excluding, during said transforming into the slant-stack domain, seismic data having a preselected characteristic. The apparatus further includes means for transforming the transformed data from the slant-stack domain into the time-space domain.
In general, in another aspect, the invention features an apparatus for anti-alias filtering seismic traces including a means for calculating a Nyquist dip pN. The apparatus further includes a means for selecting one or more dip traces having magnitudes greater than the Nyquist dip. The apparatus further includes a means for calculating an alias frequency falias. The apparatus further includes a means for applying an FFT to the selected dip traces to produce transformed selected dip traces including samples. The apparatus further includes a means for zeroing all samples in each of the transformed selected dip traces with frequencies greater than falias.
Implementations of the invention may include one or more of the following. The apparatus may include means for cosine tapering the zeroed transformed selected dip traces. The apparatus may include a three-dimensional anti-aliasing filter and the means for calculating a Nyquist dip pN may include applying the following equation:             p      N        =                                        1                          Δ              ⁢                              xe2x80x83                            ⁢                              x                2                                              +                      1                          Δ              ⁢                              xe2x80x83                            ⁢                              y                2                                                                2        ·                  f          Nyquist                      ;
where
xcex94x and xcex94y are sampling intervals in the x and y directions, respectively; and
fNyquist is the Nyquist frequency.
The apparatus may include a two-dimensional anti-aliasing filter and the means for calculating a Nyquist dip pN may include applying the following equation:             p      N        =          1                        2          ·          Δ                ⁢                  xe2x80x83                ⁢                  x          ·                      f            Nyquist                                ;
where
xcex94x is the sampling interval in the x direction; and
fNyquist is the Nyquist frequency.
The apparatus may include a three-dimensional anti-aliasing filter and the means for calculating an alias frequency falias may include applying the following equation:             f      alias        =                                        1                          Δ              ⁢                              xe2x80x83                            ⁢                              x                2                                              +                      1                          Δ              ⁢                              xe2x80x83                            ⁢                              y                2                                                                2        ·                                            p              x              2                        +                          p              y              2                                            ;
where
xcex94x and xcex94y are the sampling intervals in the x and y directions, respectively;
px is the current dip in the x-direction;
py is the current dip in the y-direction; and
falias is the Nyquist dip.
The apparatus may include a two-dimensional anti-aliasing filter and the means for calculating an alias frequency falias may include applying the following equation:       f    alias    =      1                  2        ·        Δ            ⁢              xe2x80x83            ⁢              x        ·        p            
where
xcex94x is the sampling interval in the x direction; and
p is the current dip.
In general, in another aspect, the invention features a method for p-anti-alias filtering seismic traces. The method includes applying a weighted slant stack algorithm to the seismic data for a flat dip range to produce a first result. The method further includes subtracting the first result from the seismic data to produce a second result. The method further includes applying the weighted slant stack algorithm to the second result to produce a third result. The method further includes adding the third result to the first result.
Implementations of the invention may include one or more of the following. The flat dip range may be centered around zero dip.